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New solutions of 3+ ID covariant kinetic theory are pre- 
sented for nuclear collisions in the energy domain E cm ~ 200 
AGeV. They are obtained using MPC, a new Monte-Carlo 
parton transport technique that employs very high parton 
subdivision that is necessary to preserve covariance. The 
transport results are compared with ideal hydrodynamics so- 
lutions. We show that the transport evolution differs sig- 
nificantly from hydrodynamics. In addition, we compare the 
transport freeze-out distributions to those obtained from ideal 
hydrodynamics with the Cooper-Frye isotherm freeze-out pre- 
scription. The transport freeze-out four-volume is shown to 
be sensitive to the reaction rates and deviates from both time- 
like and spacelike freeze-out 3D hypersurfaces commonly as- 
sumed. In particular, we find that there does not exist a 
universal freeze-out temperature. Finally, the transverse mo- 
mentum distributions are found to deviate by up to an order 
of magnitude from (Cooper-Frye frozen) hydrodynamics for a 
wide range of possible initial conditions and reaction rates at 
RHIC energies. 

PACS numbers: 25.75.-q, 24.10.Jv, 24.10.Lx, 25.75.Ld 
Keywords: kinetic theory, ultrarelativistic nuclear collisions, 
hydrodynamics, freeze-out, collective flow 

I. INTRODUCTION AND CONCLUSIONS 

A theoretical framework to study nonequilibrium dy- 
namics is provided by Boltzmann transport theory. The 
dynamical variables of this theory are the Lorentz- 
covariant, one-particle phase space distributions fi(x,p); 
while the dynamics is governed by transition probabili- 
ties VFc—tc' , which are Lorentz-covariant functions of the 
particle momenta. The theory, while not exact, is rather 
general. First, it is not restricted to particular parti- 
cle types. The particles could be partons, hadrons, or 
molecules. Second, the reaction rates that specify the 
dynamics are also unrestricted in their origin. For ex- 
ample, the rates could emerge from an effective quantum 
field theory or Newtonian mechanics. The primary limi- 
tations of the theory are the neglect of dynamical corre- 
lations (one-body truncation of the formal BBGKY hier- 
archy M) and the inability, without additional classical 
fields, to model phase transition dynamics. 

We consider here the simplest form of Lorentz- 
covariant Boltzmann transport theory in which the on- 
shell phase space density f(x,p), evolves with an elastic 
2^2 rate as 



P^dji = jjjihh - hh) W l2 ^ i 5\ Pl +p 2 - p 3 - Pi ) 
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+ S(x, Pl ). (1) 



Here W is the square of the scattering matrix element, 
the integrals are shorthands for J = J Aj?t5%-. , 



where 



g is the number of internal degrees of freedom, while 
fj = f(x,pj). The initial conditions are specified by the 
source function S(x, p), which we discuss later in Sec- 
tion [n|. For our applications below, we interpret f(x,p) 
as describing an ultrarelativistic massless gluon gas with 
g = 16 (8 colors, 2 helicities). 

Recall several important properties of Eq. (jl|). First, 
the particle number current and the energy-momentum 
tensor are given by 



and 



T» v (x) = 



d 3 p 



d 3 p 



jfp v f{x,p). 



(2) 



(3) 



With these definitions, particle number and energy- 
momentum conservation follow from Eq. (jl]) (when the 
source term S(x, p) = 0). Second, there is a class of fixed 
points, called global equilibria, which are phase space den- 
sities of the form 



f(x,p) 



exp 



T 



(4) 



where it M is a constant four-vector that specifies a global 
flow velocity, while T and fi correspond to the constant 
temperature and chemical potential. Furthermore, the 
H-theorem [fij states that the Boltzmann transport equa- 
tion drives the system towards global equilibrium. 

Despite its relatively simple form, the Boltzmann equa- 
tion is nonlinear with very few known analytic solutions. 
Until recently, progress to obtain even numerical solu- 
tions has been hampered by its numerical complexity. 
The rapid increase in computational power has finally 
made it possible to break through this barrier. For 
nuclear collision applications, new numerical algorithms 
are being developed, tested, and made available via the 
World Wide Web under a new Open Standard for Codes 
and Routines (OSCAR) ||. The present work is a further 
step in that development. 

For nuclear collisions at SPS energies (y/s < 20 
AGeV), numerical solutions of hadronic transport models 
have been available for some time Q . However, for higher 
collider energies, the emergence of massless partonic de- 
grees of freedom creates the technical challenge of how to 
retain Lorentz covariance. In this paper we present re- 
sults based on a new numerical technique, MPC 0.1.2 [Q, 
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that provides reliable solutions in this ultrarelativistic 
regime. 

The difficulty of the analytic treatment of the Boltz- 
mann transport equation has forced workers in the past 
to make strong simplifying assumptions. A common sim- 
plification has been to ignore the general noncquilibrium 
problem and to postulate that local equilibrium is main- 
tained at all times. In the framework of the Boltzmann 
transport theory, this choice corresponds to substitut- 
ing a local equilibrium ansatz in place of the fixed point 
global equilibria. Allowing (u„,T, fj,) to vary with the 
coordinate x M , this ansatz corresponds to 



(27T): 



exp 



H(x) -p li v/ i (x) 



T(x) 



(5) 



It is however well known that local equilibrium is not a 
solution of Eq. (Q) . The nonlinear collision term vanishes 
in local equilibrium, but p^d^f ^ in general. Only in 
the limit when all the rates go to infinity, i.e., when the 
mean free path goes to zero, can the solution approach 
local equilibrium. 

A covariant, dynamical theory can nevertheless be pos- 
tulated based on the assumption of local equilibrium. 
That is relativistic hydrodynamics |Q, which is widely 
used in heavy-ion physics ||^| to calculate observables. 
Assuming an equilibrium initial condition specified on 
a hypersurface a^ n (x), the local energy momentum and 
baryon number conservation laws 



d»T uv = 0, 



d^N E 







(6) 



reduce to (Euler) hydrodynamical equations under the 
assumptions that local chemical and thermal equilibrium 
are maintained and that dissipation (viscosity and ther- 
mal conductivity) can be neglected. In that case, the 
energy-momentum tensor and baryon current can be ex- 
pressed as T^„ = u^^e+p) - g^p and N B ^ = it M n(x) 
in terms of the local flow velocity u^ix), local pressure 
p(x), local energy density e(x), and local proper density 
n(x). The equations form a closed system if, in addition, 
the equation of state, e(p, ub), is specified. 

It is clear that the idealization of local equilibrium may 
apply, if at all, only in the interior of the reaction vol- 
ume, where the local mean free path X(x) = \/(cm{x)) 
may remain small for a while as compared to the charac- 
teristic dimensions and gradients of the system, (x) ~ 
\dfj, loge(x)| _1 . However, these assumptions are marginal 
for conditions encountered in heavy ion collisions and cer- 
tainly break down near the surface region and through- 
out the freeze-out phase. Due to longitudinal expansion 
the density decreases as 1/t until r ~ \/3R when 3- 
dimcnsional expansion rapidly increases A beyond L. 

For small departures from local equilibrium, correc- 
tions to the ideal Euler hydrodynamic evolution can be 
calculated by taking the T^ v and Ng moments of an un- 
derlying kinetic theory. To first order in X/L, the equa- 
tions reduce to the Navier- Stokes equations. The solu- 
tions to Navier-Stokes depend therefore not only on the 



equation of state, but also on the bulk and shear vis- 
cosity and thermal conductivity transport coefficients of 
the medium. While such an approach has proven use- 
ful in nonrclativistic problems and in special relativistic 
geometries, severe problems of instability and acausality 
appear when extended into the ultrarelativistic domain 

The newly formulated, covariant, parton kinetic the- 
ory technique, MPC ||, allows us to compute the highly 
dissipative evolution during the densest partonic phase 
of the reaction in a covariant manner as well as inves- 
tigate the final freeze-out dynamics. MPC is an exten- 
sion of Zhang's covariant parton cascade algorithm, ZPC 
[3f. Both MPC and ZPC have been extensively tested 
j n],[^| and compared to analytic transport solutions and 
covariant Euler and Navier-Stokes dynamics in 1+1D ge- 
ometry. A critical new element of both these algorithms 
is the parton subdivision techniqu e pr oposed by Pang 
|p|,fl2f . As shown in detail in Section III , rather high sub- 
division ~ 100 is needed to preserve Lorentz covariancc 
numerically for massless parton evolution. 

Extensions of MPC to include inelastic 2^3 partonic 
processes |l3| ] are under development, but in this paper 
we use MPC in the pure elastic parton interactions mode 
as in ZPC g. 

The aim of this work is to calculate the sensitivity of 
the evolution and freeze-out of an ultrarelativistic (mass- 
less) parton system to the transport rates and initial con- 
ditions expected in Au+Au reactions at RHIC energies 
(y/s ~ 200 AGeV). We compare the results to the (OS- 
CAR compliant) relativistic Euler hydrodynamic code 
developed by Rischke and Dumitru dJl4[]. An ideal gas 
(e = 3p) equation of state is used in that analysis. This 
work extends Refs. [jll],[l5|], focusing on the freeze-out 
problem in 3+ ID Bjorken expansion. It is less ambitious 
than for example Ref. Jl5| ] by limiting the study to mass- 
less partons. This avoids introducing yet further com- 
plications due to dissipative hadronization and hadronic 
transport effects. The initial conditions are taken from 
the HIJING multiple mini jet generator ||[l{|. 

One of the main results of this work is shown in Fig. 
|l|. To test whether ideal hydrodynamics is an adequate 
approximation of the parton transport equation (|l]), wc 
have followed the evolution both for the Boltzmann equa- 
tion (Q) and Euler hydrodynamics from the same RHIC 
initial condition from |l6| . (See Section |y| for details 
of the simulations.) Figure [l] shows the evolution of the 
transverse energy dE t / dy at midrapidity. There is a large 
difference between the transport and hydrodynamic so- 
lutions in both the 1+1 and 3+1 dimensional case, even 
for a physically extreme, 15mb, cross section. Note, that 
these 15mb curves are equivalent^] to a solution for 3 mb 



The equivalence is due to the scaling property of Eq. 
explained in Section III. 
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cross section with five times higher [[IT] initial density 
than expected with HIJING. Therefore, large deviations 
from ideal hydrodynamics are to be expected for possible 
initial conditions at RHIC. In addition, this conclusion is 
independen t of the initial system size (see explanation in 
Subsection V A ). 

This conclusion reinforces the results of Ref. | pl| , 
where it was shown that large deviations even from the 
Navier-Stokes evolution in 1+1D Bjorken expansion are 
expected for initial densities up to four times higher than 
predicted by the HIJING model @. 

A second main result of this work is shown in Fig. |^. 
We tested whether the widely used Cooper-Frye freeze- 
out prescription could "correct" final observables for the 
neglect of the early breakdown of ideal hydrodynamics. 
Figure ^ shows the final, experimentally observable p±- 
distributions divided by the thermal initial distribution. 
We find that there is an order of magnitude difference 
between transport and hydrodynamic solutions at both 
low and at high p± for realistic (~ few mb) gluonic 
cross sections. The difference is still a factor of two to 
three even for a physically extreme, 15 mb, cross section. 
In addition, increasing the radius of the initial Bjorken 
cylinder from 2 to 6 fm does not reduce the discrepancy 
between the covariant transport solutions and those of 
Cooper-Frye frozen hydrodynamics. To get closer to the 
transport pj_-distributions, one would have to choose a 
freeze-out temperature much above the commonly as- 
sumed 100 — 150 MeV range. 

The last main result of this work is illustrated in 
Fig. [| This shows that high hydrodynamic freeze-out 
temperatures that would be needed to "fit" the trans- 
port solutions are, however, inconsistent with the space- 
time freeze-out distributions of covariant transport the- 
ory. Unlike the "sharp" space-time freeze-out particle 
distributions commonly assumed using the Cooper-Frye 
freeze-out prescription, the transport theory freeze-out 
volume is four-dimensional. Particles freeze out over a 
large four-volume that forms a wedgelike freeze-out re- 
gion in the r — R plane. This freeze-out distribution de- 
pends strongly on the microscopic reaction rates (higher 
rates lead to a later freeze-out). It is not possible to tune 
the Cooper-Frye freeze-out temperature to reproduce the 
cascade freeze-out distributions. Though one can arrange 
that the Cooper-Frye freeze-out curve follows more-or- 
less the ridge, the transport distribution along that ridge 
is not correctly reproduced by Cooper-Frye frozen hydro- 
dynamics. In particular, hydrodynamic freeze-out sur- 
faces with a timelike section result in unphysical spikes 
in the dN/dr distribution which are not present in the 
transport theory calculations. These spikes arise when 
the freeze-out temperature is such that the interior of the 
system freezes o ut due to longitudinal Bjorken expansion 
(see Subsection VB 1 for further discussion). 

In summary, our results show that for a rather wide 
range of initial conditions at RHIC energies, the evolution 
of the system deviates strongly from Eulerian hydrody- 
namics throughout the 3+1D evolution. It is not possible 



to mimic the observables from the nonequilibrium evolu- 
tion by simply applying the isotherm Cooper-Frye freeze- 
out prescription to ideal hydrodynamics. The space-time 
four- volume of freeze-out, even for the largest (i? ~ 6 
fm) nuclei, does not resemble a timelike surface. In ad- 
dition, the observable transverse momentum spectra are 
very sensitive to the microscopic reaction rates. These re- 
sults indicate that while ideal hydrodynamics is a useful 
model to explore possible collective dynamics in nuclear 
collisions, the interpretation of experimental observables 
must take into account the finite transition probabilities 
W{i}_»{ji that govern the nonequilibrium evolution. For- 
tunately, numerical techniques such as MPC and ZPC are 
now readily available ||. Experimentally, the A, E cm , 
and multiplicity dependence of the observables provides 
the best way to measure these effective reaction rates in 
the ultradensc matter formed in nuclear collisions. 



II. COVARIANT PARTON TRANSPORT 
THEORY 

Equation (|l|) is the simplest form of classical Lorentz- 
covariant Boltzmann transport theory. In principle, the 
transport equation could be extended for bosons with the 
substitution /i/ 2 — + /i/2(l + /3)(l + /4) and a similar one 
for / 3 /4 (where we used the short-hand f i = f(x,p i )). In 
practice, no covariant algorithm yet exists to handle such 
nonlinearities. We therefore limit our study to quadratic 
dependence of the collision integral on /. 

The elastic gluon scattering matrix elements in dense 
parton systems are typically of the Debye-screened form: 
da/dq 2 « (9ira 2 /2)/(q 2 + y 2 ) 2 , which favors small angle 
scattering JTlj . However, the relevant transport cross sec- 
tion is a t = j da sin 2 9 cm ~ (9ira 2 /2s) log(s/4/i 2 ), where 
s w 17T 2 . In order to maximize the equilibration rate for 
a fixed cross section, we take here an isotropic differen- 
tial cross section in the center-of-mass frame instead. We 
further assume an energy-independent cross section with 
a threshold specified by /i 2 , i.e., our solutions therefore 
correspond to the microscopic dynamics specified by the 
following idealized model 



da= Q(s 2 -//)-^-dfi. 



4"7T 



(7) 



The transport cross section is 2cro/3 in this case. 

It is important to emphasize that while the cross sec- 
tion suggests a geometrical picture of action over fi- 
nite distances, we use Eq. (0) only as a convenient 
parametrization to describe the effective local transition 
probability, W. In the present study this is simply mod- 
eled as dW/dQ, — s da/dfl. The particle subdivision 
technique (see next Section) needed to recover covari- 
ance removes all notion of nonlocality in this approach, 
just like in hydrodynamics. Thus, the cross sections, e.g., 
60 mb, used in the present study to simulate rapid local 
changes of the phase space density in no way imply that 
distances bigger than 1 fm play any role. 
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With the above cutoff //, freeze-out of a test particle 
can arise in two different ways: either the system be- 
comes too dilute, i.e., l/n<7o >> L, or the system cools 
down and the threshold suppresses further interactions. 
By construction, the possibility for the latter case occurs 
along an isotherm, Tf w /x/\/l7. With Eq. (0), we can 
therefore study the influence of dissipative phenomena 
by varying the two scales o~a and /i 2 . The evolution was 
performed with uq = 3, 15, 30, 60, 121 mb and fi = 0, 
0.1, 0.5 GeV. 

The initial condition was taken to be a longitudinally 
boost invariant Bjorken cylinder in local thermal and 
chemical equilibrium at temperature T(tq) — 500 MeV 
at proper time tq = 0.1 fm/c as by fitting the gluon mini- 
jet transverse momentum spectrum predicted by HIJING 
p6| . In order to compare to hydrodynamics, we assume 
that the transverse density distribution is uniform up to 
a radius i?o = 2, 4, 6, or 8 fm. The pseudo-rapidity 
rj = l/21og((t + z)/(t — z)) distribution was taken as 
uniform between \rf\ < 5. Since we want to compare 
to chemically and thermally equilibrated hydrodynam- 
ics,^the equilibrium initial gluon density was taken for 
this T(t ) to be 



/-> f' = lf, W->W' = W/1 (a 



all). (9) 



dN 



n v ,o 



dr\d 2 x\ 



9 rrn3^ 

= — 5-J TQ 



2.65 fm" 



(8) 



Evolutions from different initial densities can be obtained 
by varying the cross section only and using the scaling 
property explained in the next Section. 



III. PARTON SUBDIVISION AND SCALING OF 
SOLUTIONS 

We utilize the parton cascade method to solve the 
Boltzmann transport equation (|l|). A critical drawback 
of all cascade algorithms is that they inevitably lead to 
numerical artifacts because they violate Lorentz covari- 
ance. This occurs because particle interactions are as- 
sumed to occur whenever the distance of closest approach 
(in the relative cm.) is d < ^cro/n, which corresponds 
to action at a distance. To recover the local character of 
equation (|l]) and hence Lorentz covariance, it is essen- 
tial to use the parton subdivision technique This 
is based on the covariance of Eq. (|l]) under the transfor- 
mation 



2 Technically, MPC was run with an out-of-chemical- 
cquilibrium initial gluon density n^o = 4 fm~ 2 as obtained 
via HIJING including final state radiation and with cross sec- 
tions (To = 2, 10, 20, 40, and 80 mb. As explained in Sec- 
tion [flj|, because of the scaling property of the solutions of 
the transport equation, the solutions for the chemical equilib- 
rium initial condition are identical when the cross section is 
rescaled by a factor I = 2.6505/4 « 1/1.509. 



As shown in Ref. p2| , the magnitude of numerical 
artifacts is governed by the diluteness of the system 
yfo I Xmfp, that scales with Lorentz violation 

therefore formally vanishes in the I — > oo limit. 



A. Convergence with subdivision 

Figure || illustrates the severeness of the cascade nu- 
merical artifacts in the case of insufficient particle sub- 
division. The top plot in Fig. ^ shows that the parton 
cascade solution for the evolution of the transverse en- 
ergy per unit rapidity does not converge until the subdi- 
vision factor reaches I ~ 100. The lack of covariance can 
be seen in the difference between the solutions in frames 
separated by 3 units of rapidity. The very fact that the 
cascade evolution is different for different particle subdi- 
visions means that the subdivision covariance (0) is itself 
violated by the cascade algorithm. Nevertheless, both 
Lorentz and subdivision covariance are recovered when I 
is sufficiently large. 

The large overshoot in the dE t jdy evolution is a result 
of the superluminal signal propagation speed inherent to 
the cascade algorithm. A cascade particle can influence 
almost instantaneously another cascade particle that is 
within the interaction range r a = \J a/ir. In a very dense 
system, a "chain" of almost instantaneous interactions 
can occur causing long range superluminal artifacts. 

As a measure of the signal propagation speed in a non- 
local collision in the cascade we define 

^ ^-partner (^collision) ^-particleiilast collision) (10) 



collision 



tlast collision 



Analytically, the deviation of the signal propagation 
speed from the speed up light can be roughly approxi- 
mated by 



Av s 



(11) 



where t is the time between the collision and the previ- 
ous collision, while r is the distance between the collid- 
ing particles at the time of the collision (r < r a ). This 
is a pessimistic estimate that maximizes the deviations. 
Assuming that subsequent collisions are uncorrelated, t 
follows a Poisson distribution^ 



P(t) 



dn 



1 



exp(— 1/\). 



dt A 

Hence the distribution of Av s (with r fixed) is 



P(Av s ) 



dn dt 



dt dAv s {Av s ) 2 X 



exp 



I A w.J A 



(12) 



(13) 



The scaling leaves the mean free path A invariant. 
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Particle subdivision reduces r CT as r a (l) = r a (V)/Vl- 
Therefore, in the large subdivision limit, the subluminal 
and superluminal tails of the signal velocity distribution 
scale as a power law 



P(Av s 



l<Ti s 



(14) 



i.e., the distribution gets narrower as Pi(Av s ) ~ 

The "measured" cascade distributions of the magni- 
tude of the signal propagation speed, dn/dv 8 , as defined 
via (|l0|), and the magnitude of its transverse component, 
dn/dv Sy ±, are shown in Fig. |J Though the distribution 
of v s is strongly peaked at v s = c, both super- and sub- 
luminal propagation are present. While increasing par- 
ticle subdivision decreases the deviations from the exact 
propagation speed c, convergence is slow. Even for a 
subdivision of 100, 13% of the collisions correspond to a 
signal propagation velocity larger than 1.5c. One must 
keep in mind that Fig. |] shows the distribution of the 
signal propagation speed measured over the length and 
time scale of a single collision. On larger scales, the de- 
viation is reduced because the large scale signal velocity 
is the sum of many small scale signal velocities. 

In summary we demonstrated that the numerical arti- 
facts due to Lorentz violation and acausality are reduced 
by subdivision and the cascade solution converges as I in- 
creases. In the I — > oo limit the cascade technique gives 
the correct numerical solution of the transport equation 
([j]). In practice, rather high subdivisions were found nec- 
essary to recover covariance. We could explore conver- 
gence up to I = 800, 200, 150, and 100, for R = 2, 4, 6, 
and 8 fm with the workstations available to us. 



B. Scaling of the transport solutions 

Subdivision covariance (|9|) actually implies that the 
transport equation has a broad dynamical range, and 
the solution for any given initial condition and transport 
property immediately provides the solution to a broad 
band of suitably scaled initial conditions and transport 
properties. This is because solutions for problems with 
/ times larger the initial density dN/drjd 2 x±, but with 
one l-th the reaction rate can be mapped to the original 
(1 = 1) case for any I. We must use subdivision to elimi- 
nate numerical artifacts. However, once that is achieved, 
we have actually found the solution to a whole class of 
suitably rescaled problems. 

The dynamical range of the transport equation (|l|) 
is further increased by its covariance under coordinate 
rescaling 



f(x,p) -» f(x,p) = f[ -,p 



w 



(15) 



This is a simultaneous rescaling of space-time and the 
transition probability. In addition, there is also a covari- 
ance under rescaling of the momenta 



w({pi})^w'Upi}) = i*w 



(16) 



such that the particle density is again unchanged. This 
scaling also implies a rescaling of the mass m — > ml — 
m/lp. Combining the three scaling transformations, we 
find covariance of the transport theory under 



f(x,p)^ f'(x,p) = l-Hf 



x p 

F' F 



w({ Pi }) 



(17) 



In our calculation using MPC, we vary the physical 
parameters: a, fi, Tq, Rq, tq, and n v ^ = dN / dnd 2 x±_\ To 
(the rapidi ty i nterval r/max = 5 was fixed). Keeping in 
mind Eq. (A3) and that 



dN 



dyd% 



J d 2 p ± drj TO t ch(y-r/)r/(x_L , T), t, p±,y), 



covariance under the transformation (|T^) implies that 
once the solution for a particular choice of these param- 
eters is known, then the solution is known for any other 
choice of the parameters which are related to the original 
via 



a ' — lx ^ 1<T 7 Pq ~ ^p^Oj R'q — IxRoi 



n' vQ = l x ln n . Q , /j,' = Ipfi, 



IxTQ. 



(18) 



Therefore, we can scale one solution to others provided 
that h/Tq, Rq/tq, and ero, ; o ~ tq/Xmfp remain the 
same. For example, three times the density with one- 
third the cross section leaves all three parameters the 
same, hence the results can be obtained via scaling with- 
out further computation. Table [j] shows sets of the three 
ratios that we mapped out via MPC. 



IV. FREEZE-OUT 

A. Hydrodynamic Freeze-Out Problem 

In Section | we argued that hydrodynamics cannot be 
valid during the complete evolution in nuclear collisions 
because the assumption of local equilibrium breaks down. 
Thus, in spite of its appeal, hydrodynamics cannot be 
compared with measurements without additional model 
assumptions needed to specify when and how it breaks 
down. The problem of determining those extra model 
assumptions is the so-called freeze-out problem. 
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For application to nuclear collisions, freeze-out cannot 
be formulated as an expansion in X/L since by definition 
it occurs when that ratio exceeds unity. Hence, even the 
Navier-Stokes hydrodynamics is inadequate to solve the 
freeze-out problem. 

A common freeze-out prescription, which we here name 
"Cooper-Frye frozen hydrodynamics" , is to assume the 
validity of ideal hydrodynamics up to a "sharp" 3D 
freeze-out hypersurface a^(x). Assuming that all inter- 
actions suddenly cease on that hypersurface, the final 
(frozen-out) invariant differential distribution of particles 
is then computed via the Cooper-Frye formula [[L8l : 



EdN = 



d 3 p 



da M (x)p li f(x,p). 



(19) 



Here dcr M (x) is the normal vector to the 3D freeze-out 
hypersurface at the point x, while f(x,p) is assumed to 
be in local equilibrium and hence, for classical particles, 
given by Eq. (|^). While this prescription is covariant 
and appealingly simple, it suffers from several well known 
problems 

First, because the hydrodynamical solutions do not 
contain dynamical information needed to compute the 
freeze-out hypersurface, the assumed one is simply an ad 
hoc external constraint. It is usually parameterized in 
terms of a few physically "reasonable" parameters, the 
most common being a freeze-out isotherm T(a fJ ') = Tj 
or freeze-out energy density e/ . It is not possible to esti- 
mate the errors introduced by such a prescription. 

Second, the Cooper-Frye formula allows negative con- 
tributions to the measurable particle yields ^ |22| . This 
can be avoided pl|-p3| by choosing a non-equilibrium 
post freeze-out distribution that does not have particles 
in the phase space domain where da^p^ < 0. However, 
such a choice still relies on the existence of a sharp 3D 
freeze-out surface. 

Finally, while an idealized sharp freeze-out surface may 
be adequate for applications to quasi-stationary macro- 
scopic systems, it cannot be justified in expanding meso- 
scopic systems in which L/Xis never large. The very fact 
that such systems do freeze out, i.e., A(cr) > L, means 
that the solution to freeze-out problem must entail global 
information as the system becomes more and more di- 
lute. Furthermore, there is no way to justify the neglect 
of final state interactions during freeze-out stage of the 
reactions while expansion and rarefaction are causing the 
system to depart from local equilibrium. 

In Ref. a continuous emission hydrodynamical 
freeze-out model was proposed to overcome some of these 
problems. The global information relevant to freeze-out 
in that model is taken there as the Glauber escape prob- 
ability 

P(x,p) = exp [ - / dr' av re in(x{T)) ] . (20) 



x^, the line integral runs over the future trajectory and 
therefore is exponentially sensitive to the future evolution 
of the system. This leads to a formidable self-consistency 
problem. For special geometries such as Bjorkcn boost 
invariant expansion, a rough estimate of P can be made 
using the approximate scaling Bjorken hydrodynamic so- 
lution 

n{x{r)) = ^n(x(r )) Q(R 2 - (x ± + v ± r) 2 ). (21) 

T 

Together with the Glauber straight line trajectory, this 
leads to the characteristic power law survival probability 



P~(-L. 

\Tout 



UT n(T ) 



(22) 



This formula reveals clearly the highly nonlocal charac- 
ter of the freeze-out problem. At any point in spacetime, 



that also appears, ej* in the J/ip suppression problem 
p5| . While Eq. (^(]) captures essential global physics 
of freeze-out, it is not complete since before freeze-out 
the trajectories cannot be straight if local equilibrium 
is maintained via the assumed hydrodynamic equations. 
Also, in the surface region where P ~ 1/2 neither hydro- 
dynamics nor eikonal dynamics applies. 

The solution to the freeze-out problem in classical me- 
chanics is given by microscopic transport theory. A hy- 
brid approach that partially reaches that end was pro- 
posed in Ref. [ p6[ , which combines partonic hydrodynam- 
ics with hadronic transport theory. In that approach, 
hydrodynamics is assumed to hold up to only some in- 
termediate critical temperature hypersurface, T{a^ nt ) = 
T c > Tf, on which the fluid is converted to hadrons 
via the Cooper-Frye formula. Subsequent evolution of 
the hadronic system is then calculated by solving the 
hadronic transport theory as encoded in UrQMD. As 
noted in Ref. |24|, the freeze-out surface is actually a 
diffuse four-volume, and in addition different hadronic 
species freeze out over different four-volume domains. 
This sequential freeze-out leads to strong observable cor- 
relations such as the mass dependence of final transverse 
spectra. The main limitation (and/or advantage) of the 
above hybrid model J27j is the need to assume the va- 
lidity of hydrodynamics in the dense partonic phase of 
the collision. It is advantageous in that possible collec- 
tive effects due to the quark-gluon confinement transi- 
tion can be explored with hydrodynamics using "realis- 
tic" equations of state. It is disadvantageous in that it 
is far from clear that local equilibrium is ever reached 
during the evolution. Recall [|llj that dissipative effects 
on even global observables such as the transverse energy 
per unit rapidity cannot be accurately calculated using 
the Navier-Stokes equations. 

Despite these known complications of the freeze-out 
problem, ideal hydrodynamics and Cooper-Frye freeze- 
out are still commonly used to fit experimental data 
using isotherm freeze-out hypersurfaces and draw infer- 
ences about the underlying dynamics. The consistency 
and significance of interpretations based on such fits can 
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only be assessed by comparing detailed dynamical trans- 
port calculations to the hydrodynamic limit (see Section 

0)- 

B. Formal Definition for Freeze-out 



While neither normalized nor unique, this expression pro- 
vides at least a formal definition of the freeze-out distri- 
bution for the Boltzmann equation solely in terms of the 
phase space distribution f(x,p), and the assumed tran- 
sition probabilities Wij^ki ■ 



An important experimental observable aspect of the 
space-time evolution of kinetic theory is the freeze-out 
distribution. In the framework of discrete parton cas- 
cade dynamics, the definition of the freeze-out distribu- 
tion, dNf , is the number of partons per d A xd A p invari- 
ant phase space volume with momentum p^ that have a 
collision at x 11 but suffer no more collisions. Given the 
trajectories, (x a (t), p a (f)) or the world lines x^{r) of all 
partons a, that distribution is given by the ensemble av- 
erage of the space-time coordinates, x^ f = (t a f, x a (t a /)), 
of the last collision together with the final outgoing mo- 
mentum, p(t a f + + ): 



ffo(x,p) 



dN 



fo 



d 4 xd 4 p 



J2S(t-t af )S 3 (x-x a (f a/ )) 



x 5 i (p~ Pa {t af + 0+))) . (23) 



Because the freeze-out times, t a f, are distributed over a 
broad time interval, ff does not correspond to 

f(x,p)=NQ dT5\x-x(r))8\p^p{T))^ (24) 

at any time or on any fixed 3-D hypersurface. Note that 
/ measures the phase space density of the world lines 
x^(t) and their four- velocities at a single point x^ . On 
the other hand, ff measures the phase space density 
of last scattering events, where the momentum p of a 
particle was last changed. Pion interferometry |2£| mea- 
sures the Fourier transform of // G . Note that even after 
integrating over the freeze-out points, the final observed 
momentum spectrum, J d 4 xff (x,p), is only equal to the 
Cooper-Frye formula if x^f happen to lie on a sharp 3D 
hypersurface, a 11 (Ci, C27 (3)- We can write the Cooper- 
Frye freeze-out distribution then as 



dN CF 
d 4 xd 4 p 



N d\6 4 (x - a(())5 4 (p - p(a(C))). (25) 



As discussed in the Appendix, we can write the (on- 
shell) freeze-out distribution in terms of the solution of 
the Boltzmann equation as 



dF fo {x,p 1 ) 
d 4 xd 3 pi 



Po(x,pi 



W34^15<5 4 (P3 +Pi - pi -Ps)hf4 
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(26) 



V. NUMERICAL RESULTS 
A. Kinetic versus Hydrodynamic Evolution 

To test the ideal hydrodynamical assumptions against 
transport theory, it is essential to eliminate as many 
model differences as possible. For example, both the hy- 
drodynamic model and the kinetic theory should have 
the same degrees of freedom, the same equation of state, 
and the same initial conditions. Equation (Q) describes a 
gas that in thermal equilibrium has the equation of state 
e = 3p, if the partons are massless. We therefore used 
this ideal gas equation of state in the hydrodynamical 
simulations. We also chose the transport initial condi- 
tion to be in local equilibrium, since hydrodynamics is 
limited to such initial conditions. 

The hydrodynamic algorithm used [fl4|] is furthermore 
designed for particles without a conserved charge, i.e., the 
particle number changes as dictated by chemical equilib- 
rium. The algorithm solves the energy-momentum con- 
servation equation to obtain the energy density, pressure 
and flow evolution. Then, instead of the charge con- 
servation equation, it exploits the relation between den- 
sity, particle mass, and temperature in chemical equilib- 
rium to compute the freeze-out particle distribution. It 
is important to note therefore that Eq. (|l|) with elas- 
tic collisions has the same hydrodynamic limit as the 
hydrodynamic model only if the partons are massless. 
This is because ideal hydrodynamics conserves entropy 
p9[ and for massless particles in thermal and chemical 
equilibrium entropy conservation is equivalent to particle 
number conservation.^ For massive particles, we would 
have to compare transport to hydrodynamics with par- 
ticle conservation. Conversely, we would need to sup- 
plement Eq. (fy) to include inelastic channels, such as 
2 <-> 3 in Ref . |13|| , to compare to chemically equilibrated 
hydrodynamics. In the infinite rate limit we recover the 
hydrodynamic model even though we have a fixed num- 
ber of particles. However, when the solution is out of 
equilibrium (either thermal, chemical, or both), it does 
make a difference whether we include particle number 
changing processes or not. 

To test whether ideal hydrodynamics is an adequate 
description of the parton transport theory ([!]) , we com- 
pare the evolution of the transverse energy dE t /dy at 
midrapidity from the two models. This comparison is 



Because in this case s = 4n. 
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free from any hydrodynamic freeze-out prescription be- 
cause the transverse energy is given directly by the phase 
space distribution as 



dE t 
dy 



d 2 p±dr] d 2 x±m t cosh(y — rf) 

x m t f(y,p±,r),x±,T), (27) 



where, through the local equilibrium ansatz ([5]), the hy- 
drodynamic phase space evolution is determined by the 
evolution of the flow velocity and local temperature as 
dictated by the equations of motion (||). 

Figure [H shows the transverse energy evolution from 
transport theory and hydrodynamics, for an initial 
Bjorken cylinder radius of 2 fm, with To = 0.1 fm/c, 
Tq= fi = 0.5 GeV, n v , = 2.6505 fm~ 2 (via scaling), a = 
15, and 60 mb. (We chose rj m = 5, subdivisions 800 for 
3+1D, 256 for 1+1D, and a 100 fm 2 transverse area for 
the 1+1D evolution.) 

The transverse energy decreases much faster from ideal 
hydrodynamics than from kinetic theory, both in 1+1D 
and 3+1D, showing that hydrodynamics does more work 
than the cascade. This is due to the different phase space 
evolution in the two models. The early discrepancy, even 
for cross sections as extreme as 15 or 60mb, indicates that 
either the transport evolution gets very quickly out of 
equilibrium, or the initial evolution is close to equilibrium 
but the energy-momentum tensor is not ideal. Note that 
even if the latter is true, it does not necessarily mean 
that this initial, locally equilibrated, nonideal dynamics 
can be described by the Navier-Stokes equations. 

The above conclusion holds for any initial system size 
larger than 2fm as well. Since the 1+1D curves corre- 
spond to the infinite transverse size limit, the hydro- 
dynamic and transport evolutions for initial sizes larger 
than 2fm will lie between the 2fm and the 1+1D curves 
for hydrodynamics and for transport theory, respectively. 
Because these two regions do not overlap, the discrepancy 
between ideal hydrodynamics and transport theory will 
not disappear with increasing system size. 



B. Kinetic vs Hydrodynamic Freeze-out Results 

In the previous Subsection we showed that parton ki- 
netic theory does not reduce to ideal hydrodynamics for 
initial conditions at RHIC. Thus, the final observables 
from the two models can be similar only if the hydrody- 
namic freeze-out prescription helps mimic the observables 
from the nonequilibrium transport evolution. 

Here we test whether one can reproduce the transport 
observables by a suitable choice of the hydrodynamic 
freeze-out parameters. We chose the widely-used Cooper- 
Frye freeze-out prescription (|l^) with isotherm freeze-out 
surfaces, despite all known problems discussed in Section 
[Tv| . Hence, our only adjustable parameter is the freeze- 
out temperature. Since Eq. (n[) describes Boltzmann 



classical particles, we must use the classical distribution 
(|J) in the Cooper-Frye formula. 



1. Coordinate space evolution 

Freeze-out distributions in space-time from MPC are 
shown in Figs. || and ||. Due to the assumed cylindrical 
symmetry and longitudinal boost invariance, that distri- 
bution is only a function of r and R. 

Figures ||, |B| show the freeze-out distribution for initial 
radii 6 fm and 2 fm, respectively, with To = 0.1 fm/c, 
T = (i = 0.5 GeV, n V! o = 2.6505 fm~ 2 (via scaling), a = 
3, 15, and 60 mb. For comparison, three different freeze- 
out isotherms are also shown from solution of Cooper- 
Frye frozen ideal hydrodynamics. (We chose r/ m = 5, 
subdivisions 800 for 2 fm, and 150 for 6 fm.) 

Unlike the sharp hydrodynamic freeze-out surface, the 
freeze-out distribution from the cascade is a broad wedge. 
Particles originate from a hypervolume in space-time, 
rather than from a hypersurface. In the top left plot 
(3 mb, 6 fm) in Fig. the wedge moved down to r = To, 
which is a general feature for very low reaction rates. In 
the limit of a vanishing reaction rate, all particles freeze 
out from t = To. 

Figures ||, || show that particles freeze out later with 
increasing microscopic rates as expected. The maximum 
of the wedge moves outward with increasing rates, hence 
no freeze-out temperature can be universal. If we tune 
the freeze-out temperature to get as close as possible to 
the cascade freeze-out distribution, the freeze-out tem- 
perature will depend on the reaction rate. 

Thus, the remarkable agreement seen in the bottom 
figure of Fig. ^ between Cooper-Frye frozen ideal hydro- 
dynamics with a 130 MeV freeze-out temperature and 
the cascade with a — 60 mb is a mere coincidence; higher 
rates would lead to a later freeze-out. For very high re- 
action rates, the 130-MeV hypersurface from the cascade 
would be very close to that from hydrodynamics because 
the hydrodynamic evolution is the infinite reaction rate 
limit of the cascade evolution. But that does not mean 
that the freeze-out distributions are the same. On the 
contrary if hydrodynamics and the cascade are close to 
each other at T = 130 MeV then we have no justification 
to stop the hydrodynamic evolution and freeze out with 
Eq. ( |H| ) because we are still in equilibrium and particles 
will certainly collide in the future, i.e., they have not yet 
frozen out. 

It is not possible to tune the Cooper-Frye freeze-out 
temperature to reproduce the cascade freeze-out distri- 
bution. Though the contour plots in Figs. 0, || suggest 
that such a tuning can get the hydrodynamic freeze-out 
curve close to the ridge of the wedge of the cascade freeze- 
out distribution, that is not enough. As the dN/dr dis- 
tributions show, the resulting hydrodynamic freeze-out 
distribution is not close to the cascade distribution be- 
cause one has to reproduce not only the curve given by 



8 



the ridge of the wedge but also the exact distribution 
along this curve. 

If the freeze-out temperature is high enough to yield 
a freeze-out surface with a timelike portion, we get un- 
physical spikes in the freeze-out distribution that are 
not present in the cascade calculations. This can be 
seen in Fig. § for T f = 200 MeV, and in Fig. | for 
T f = 130 and 200 MeV. For example, for R = 6 fm 
with Tf = 130 MeV, Cooper-Frye frozen hydrodynam- 
ics produces most particles at around r = 5.6 fm/c. 
Cooper-Frye frozen hydrodynamics produces most par- 
ticles at around r = 5.6 fm/c. This is because the in- 
side of the cylinder follows a ID Bjorken evolution with 
T(t) — ^(ro/r) 1 / 3 until the rarefaction wave from the 
boundary arrives. The rarefaction wave travels with a 
speed c s — 1/V3. If the system is large enough, most of 
the system reaches the freeze-out temperature before the 
rarefaction wave arrives, i.e., during the 1+1D Bjorken 
evolution. With our parameters T — 0.5 GeV, Tf Q = 130 
MeV, and To = 0.1 fm/c, this gives a freeze-out for the 
inside of the cylinder at Tf = 5.6 fm/c, which is in com- 
plete disagreement with our transport theory solutions. 
Furthermore, it does not correspond to the infinite reac- 
tion rate limit either because in that case particles freeze 
out very late. 

Hence the peaks in dN/dr at r = 5.6 fm/c (Tf — 
130 MeV) and r = 1.6 fm/c (Tf = 200 MeV) are a 
clear consequence of the arbitrary freeze-out prescription 
using Eq. ( |l9| ) with isotherm freeze-out hypersurfaces. 
Smearing the peaks out around their maxima does not 
help either because that does not change the location of 
the peaks, while the maximum from the cascade moves 
outward with increasing reaction rates. 

2. Momentum space 

The freeze-out distribution in momentum space is 
shown in Figs. || and [| Figure [| shows the freeze-out 
p_i_-distribution for initial radii 2 fm and 6 fm, cascade 
cross sections 3, 15, and 60 mb compared to ideal hydro- 
dynamics with a Cooper-Frye freeze-out at temperatures 
Tf = 100, 130, and 200 MeV. As the reaction rate in- 
creases, the small pj_-slopes rise as the system cools due 
to longitudinal work. The p^-distribution seems to ap- 
proach that of Cooper-Frye frozen hydrodynamics. How- 
ever, this is only an illusion on a low-resolution logarith- 
mic plot. Figure]^, where we plotted the final pj_-spectra 
divided by the initial Tq — 500 MeV thermal one, shows 
that there is a large, up to a factor of ten difference at 
both low (< 0.5 GeV) and high p x (> 2 GeV), depend- 
ing on the microscopic rates. For all the cases studied, 
Cooper-Frye frozen hydrodynamics has more low-pj_ par- 
ticles but fewer high-p^ ones than the cascade. This is 
not necessarily a general feature because the assumed hy- 
drodynamic freeze-out temperature is an arbitrary num- 
ber. A later freeze-out (lower temperature) gives a larger 



slope, an earlier freeze-out (higher temperature) gives a 
smaller one. 

It is also striking that one would need rather high, 
Tf ~ 300 — 450 MeV freeze-out temperatures to get closer 
to the cascade pj_-spectra. We conclude that it is not pos- 
sible to reproduce both the space-time and the momen- 
tum space transport theory freeze-out distributions us- 
ing ideal hydrodynamics with the isotherm Cooper-Frye 
freeze-out prescription. Either one needs to treat hy- 
drodynamic freeze-out more accurately than the Cooper- 
Frye prescription, or one needs to use full-scale transport 
theory instead of ideal hydrodynamics. The present work 
is a step in the latter direction, while Refs. [21-23] are 
important steps in the former direction looking for a sim- 
plification of the full transport theoretical problem that, 
hopefully, will still be applicable to a wide class of situa- 
tions. 



VI. OUTLOOK 

There are many open problems in the development of 
covariant transport theory. The most urgent need is to 
develop practical convergent algorithms to incorporate 
inelastic 2^3 processes to allow studies of chemical 
equilibration. Preliminary work in Ref. |]l3| indicated a 
rather slow convergence towards Lorentz covariance with 
particle subdivision. Unlike the Z -1 / 2 convergence in 2 — > 
2, a much slower ~ Z" 1 / 5 convergence is expected in 2 <-> 
3 processes even when nonlocal formation physics (At > 
h/AE) is neglected. 

Also, we note that all results in this paper pertain 
to homogeneous initial conditions. In Ref. p0| , it was 
shown that jets induce large nonstatistical local fluctua- 
tions that may evolve in a turbulent manner. A transport 
study of the evolution from such inhomogeneous initial 
conditions would be useful to compare to the known hy- 
drodynamic solutions. 
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APPENDIX: FORMAL DEFINITION FOR 
FREEZE-OUT 

Unlike in the cascade solution where the freeze-out dis- 
tribution is trivially defined by Eq. (|25|), in the Boltz- 
mann equation / changes continuously and no discrete 
final collisions can be identified. In this appendix we pro- 
pose a generalization of Eq. (|20| ) which is independent 
of the discrete numerical cascade picture. We motivate 
here a formal definition, Eq. ( ptf ) , of the freeze-out dis- 
tribution using solely f(x,p) and Wij^kl- 

Following the notion of the "last collision", one can 
first compute the probability that a particle starting at 
a coordinate Xi with momentum pi does not collide any 
further. The collision rate is given by 

-(x,pi,p 2 ) 



cou - d±x 

= fi(x,Pi)f2(x,p 2 )o-(pi,P2)vi2d 3 p 1 d 3 p2, (Al) 

where the relative velocity and the total cross section are 
given with the Lorentz scalar 



tu = \J {PiP2^) 2 - m\m\ 



(A2) 



as 



V\2 



tl2 

E\E 2 ' 



v{PI,Pi) = — JJw i2 ^34S 4 {Pl +P2 ~P3 -Pi)- (A3) 
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A free pointlike particle has the phase space distribu- 
tion 



fi(x,p)= drS' 1 (x-x 1 -u 1 T)S' 1 (p-pi), (A4) 
Jo 

where for an on-shell particle 

S (p° - \Jp 2 + m 2 ) = 5(p 2 - m 2 )2\/m 2 + p 2 @(p ), 
i.e., 



8\p - Pl ) = 5 3 (p - Pl )S( P 2 - m 2 )2^firf + pje(p ) 
and thus^J 



/i(x,p)= / drd (x-Xi-Uxt) 
Jo 



X 3, \ V™ 2 + Pl 

x S J (p- pi)-^ L -. 



(A6) 



Plugging this result into Eq. (Al), the probability that 
a free particle will not have any further collisions isQ 



Po(zi,Pi) = exp 



= exp 



drd 3 p2 
E 2 m 



T co ud 3 pid 3 p2d 4 x 



f2(xi + UlT,p2)o-(p 1 ,p 2 )tl2 ■ (A7) 



Now we can write the freeze-out distribution as the 
number of particles having a collision at x^ with outgoing 
momentum pi times the probability that these particles 
do not collide any further, i.e., 



Ei 



dFff(x, Pl ) 



-Po(x,pi) 



d 4 xd 3 pi 

x 2 ///w 3 4^i5<5 4 b 3 +P4-Pi-p 5 )/3/4. (A8) 
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This definition does not include those particles that 
are formed but suffer no collisions afterward. Their con- 
tribution is 



Ey 



dF 



form 



Jo 



(x,Pi) 



d A xd?p\ 



S(x,pi)P (a;,Pi)- 



(A9) 



Hence the final freeze-out distribution is given by Eq. 

The definition ( |26| ) should be regarded only one mea- 
sure of the freeze-out distribution because it has several a 
shortcomings. The probabilities summed are not proba- 
bilities for disjoint events. One should exclude the volume 
in space time given by all the linear paths of the already 
frozen-out particles. This requires knowledge of multi- 
particle correlations beyond the scope of the Bolztmann 
equation. As long as those excluded volume effects are 
small, Eq. ( |26| ) is adequate. A clear problem with the 
present formal definition is that particle number and mo- 
mentum are not conserved by it as is automatic in (|2^). 
It is interesting to contrast on the other hand, the triv- 
ial way that the cascade solves this problem through Eq. 
(p3|). In cascade, the TV-body correlations are automati- 
cally calculated and freeze-out is easily defined conserv- 
ing number and total four-momentum. The continuum 
limit is thus subtle. Our numerical results define that 
continuum limit as the limit of infinite subdivisions us- 
ing the cascade technique. 



5 Recall, the on-shell phase space distribution is defined via 

f(x,p) = 2m8(p 2 - m 2 )e(p )/(x,p). (A5) " This can be derived assumin S that subsequent collisions 

are uncorrelated (just like the similar formula for the inho- 
mogeneous Poisson distribution). 
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MPC vs hydro (1+1D and 3+1D) 
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FIG. 1. This figure shows the evolution of the transverse en- 
ergy dEt/dy at midrapidity, normalized by the initial transverse 
area, from kinetic theory and from hydrodynamics both for 1+1 
(transverse periodic) and 3+1 dimensions. The initial distribution 
was a Bjorken cylinder with a radius Ro = 2 fm at proper time 
tq = 0.1 fm/c in local thermal and chemical equilibrum at To = 500 
MeV. The cross sections were a = 15 and 60 mb, with the cutoff 
fi = 0.5 GeV. Note, that the hydrodynamical results are free from 
any arbitrary freeze-out prescription because they depend only on 
the evolution of the phase space distribution ansatz (^) as dictated 
by the equations of motion (^j). The difference between the hydro- 
dynamical and the kinetic theory results, even for a large, 60 mb 
cross section, indicates that ideal hydrodynamics is not applicable 
for the evolution. 
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FIG. 2. The top figure shows the freeze-out px-distributions 
relative to the initial (To = 500 MeV) p± -distribution from an ini- 
tial radius Ro = 2 fm for the cascade with a = 3, 15, 60 mb, and for 
ideal hydrodynamics with Cooper-Frye freeze-out with freeze-out 
temperatures T f =100 MeV (dashed-dotted line), 130 MeV (thick 
solid line), and 200 MeV (dotted line). The bottom figure shows 
the same but for an initial radius Ro = 6 fm with cross sections 3 
and 15 mb. 
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FIG. 3. The left column shows the transverse coordinate and 
proper time distribution R ^ dT , of frcczc-out coordinates. Top 
row corresponds to a = 3 mb and middle row to 15 mb. The 
initial Bjorkcn cylinder radius is Rq = 6 fm in both cases. The 
right column shows contour plots corresponding to the left col- 
umn. The thick lines show Coopcr-Fryc isotherms: Tf = 100 MeV 
(dashed-dottcd line), 130 McV (thick solid line), and 200 MeV (dot- 
ted line). The bottom figure compares the proper time freeze-out 
distribution, 4^, for the different cases. 
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FIG. 4. The top figure shows the Et evolution at midrapidity 
(solid lines) for subdivision factors 2 = 1, 10, 100 from the initial 
condition Rq = 2 fm, a = 10 mb, n V} o = 4 fm~ 2 , To = £t = 0.5 
GeV, to = 0.1 fm/c. The dashed curves are for the same initial 
condition but all particles were boosted longitudinally by 3 units of 
rapidity and dEt/dy was computed at y = 3. Though not visible, 
the curve for 2 = 400 falls on top of the 2 = 100 curve. The bottom 
plots show the normalized distribution of the signal propagation 
velocity averaged over the full evolution. 
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FIG. 5. The left column shows the transverse coordinate and 
proper time distribution R ^ dT , of freeze-out coordinates. Top row 
corresponds to a = 15 mb and middle row to 60mb. The initial 
Bjorkcn cylinder radius is Rq = 2 fm in both cases in contrast to 
Fig. js| where Rg = 6 fm. The right column shows contour plots 
corresponding to the left column. The thick lines show Cooper-Frye 
isotherms: T f = 100 MeV (dashed-dotted line), 130 MeV (thick 
solid line), and 200 MeV (dotted line). The bottom figure compares 
the proper time freeze-out distribution, s££, for the different cases. 
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FIG. 6. The top figure shows the freeze-out p± -distributions 
from an initial radius Ro = 2 fm for the cascade with a = 3, 15, 
60 mb, and for ideal hydrodynamics with Cooper-Frye freeze-out 
with freeze-out temperatures Te =100 MeV (dashed-dotted line), 
130 MeV (thick solid line), and 200 MeV (dotted line). The initial 
px -distribution is shown using pluses. The bottom figure shows 
the same but for an initial radius Ro = 6 fm with cross sections 3 
and 15 mb. 
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TABLE I. Solutions for the above sets of parameters were 
computed via MPC for the present study. 



